% x: sorting variable (size, value); y: sorted variable
% (return);u:frequency of sorting; m: number of portfolios
% number: how many elements in each portfolio
% belong: which portfolio each element belongs to

% the length of x is the length of y divided by u
function[st,number,belong,resi]=sorting(x,y,u,m)
% set space
n=length(y(1,:)); % size of cross section
T=length(y(:,1)); % length of time series
st=zeros(T,m);% st is the sorted portfolio returns, with the same length as y
z=floor(n/m)+1; % at most, how many assets in each portfolio
% start
freq=fix(T/u); % the frequency of rebalancing
resi=T-freq*u;
if resi==0
else
    y=y(1+resi:end,:); % adjust the sample such that rebalancing is exact
end
% set space for record
for i=1:m
belong{i}=NaN*ones(freq,z); 
number{i}=NaN*ones(freq,1);
end

for i=1:freq
     d=isnan(y(u*(i-1)+1,:));% the first obs on y, for example, return for Jan
     x(i,find(d==1))=NaN; % get rid of those countries' x without y (first obs in the u period,
     % from (i-1)*u to i*u) data
     [a,b]=sort(x(i,:)); % a is the sorted data, b is the order
     no2=1-isnan(a); no=sum(no2); %no:the number of non-NaNs;
     l=NaN*ones(m,1);
     j=1; % starting from portfolio 1
     mm=m; % how many portfolios needs to be allocated
     while j<m+1 % stop when j is about to become m+1
     l(j,1)=floor(no/mm); % how many assets in portfolio j
     no=no-l(j,1); % after sorting into portfolio j, how many assets are left
     mm=mm-1;
     j=j+1; % next portfolio
     end


     
     b2=NaN*ones(m,z);
     l2=NaN*ones(m,1);
     j=1;
     while j<m+1
         if j==1
         l2(j,1)=l(1,1); 
         b2(j,1:l(1,1))=b(1:l(1,1));
         else
         l2(j,1)=l2(j-1,1)+l(j,1);
         b2(j,1:l(j,1))=b(l2(j-1,1)+1:l2(j,1));
         end
     y2=y(u*(i-1)+1:u*i,b2(j,1:l(j,1)));
     belong{j}(i,1:l(j,1))=b2(j,1:l(j,1));
     st(u*(i-1)+1:u*i,j)=nanmean(y2,2);
     c=1-isnan(y2); % how many assets in each portfolio
     number{j}(i,1)=sum(c(1,:));
     j=j+1;
     end
end

end